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Abstract 

We summarise research reported in (Efroimsky 2002, 2003; Efroimsky & Goldreich 



2003a, b) and develop its application to planetary equations in non-inertial frames. 
■ Whenever a standard system of six planetary equations (in the Lagrange, Delaunay, 

or other form) is employed, the trajectory resides on a 9(N-l)-dimensional submanifold 
of the 12(N-l)-dimensional space spanned by the orbital elements and their time deriva- 
tives. The freedom in choosing this submanifold reveals an internal symmetry inherent 
in the description of the trajectory by orbital elements. This freedom is analogous to 
the gauge invariance of electrodynamics. In traditional derivations of the planetary 
equations this freedom is removed by hand through the introduction of the Lagrange 
constraint, either explicitly (in the variation-of-parameters method) or implicitly (in 
the Hamilton- Jacobi approach). This constraint imposes the condition that the orbital 
elements osculate the trajectory, i.e., that both the instantaneous position and velocity 
be fit by a Keplerian ellipse (or hyperbola). Imposition of any supplementary con- 
straint different from that of Lagrange (but compatible with the equations of motion) 
would alter the mathematical form of the planetary equations without affecting the 
physical trajectory. 

For coordinate-dependent perturbations, any gauge different from that of Lagrange 
makes the Delaunay system non-canonical. In a more general case of disturbances 
dependent also upon velocities, there exists a "generalised Lagrange gauge" wherein 
the Delaunay system is symplectic (and the orbital elements are osculating in the phase 
space). This gauge reduces to the regular Lagrange gauge for perturbations that are 
velocity- independent . 

We provide a practical example illustrating how the gauge formalism considerably 
simplifies the calculation of satellite motion about an oblate precessing planet. 
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1 Introduction 



1.1 Prefatory notes 

On the 6-th of November 1766 young geometer Joseph-Louis Lagrange, 1 invited from Turin 
at d'Alembert's recommendation by King Friedrich the Second, succeeded Euler as the Di- 
rector of Mathematics at the Berlin Academy. Lagrange held the position for 20 years, and 
this fruitful period of his life was marked by an avalanche of excellent results, and by three 
honourable prizes received by him from the Academie des Sciences of Paris. All three prizes 
(one of which he shared with Euler) were awarded to Lagrange for his contributions to celes- 
tial mechanics. Among these contributions was a method initially developed by Lagrange for 
his studies of planet-perturbed cometary orbits and only later applied to planetary motion 
(Lagrange 1778, 1783, 1788, 1808, 1809, 1810). The method is based on an elegant mathe- 
matical novelty, the variation of parameters emerging in solutions of differential equations. 
The first known instances of this tool being employed are presented in a paper on Jupiter's 
and Saturn's mutual disturbances, submitted by Euler to a competition held by the French 
Academy of Sciences (Euler 1748), and in the treatise on the Lunar motion, published by 
Euler in 1753 in St. Petersburg (Euler 1753). However it was Lagrange who revealed the full 
power of the method. 

Below we shall demonstrate that the equations for the instantaneous orbital elements 
possess a hidden symmetry similar to the gauge symmetry of electrodynamics. Derivation of 
the Lagrange system involves a mathematical operation equivalent to the choice of a specific 
gauge. As a result, trajectories get constrained to some 9-dimensional submanifold in the 12- 
dimensional space constituted by the orbital elements and their time derivatives. However, 
the choice of this submanifold is essentially ambiguous, and this ambiguity gives birth to an 
internal symmetry. The symmetry is absent in the 2-body case, but comes into being in the 
N-body setting (N > 3) where each body follows an ellipse of varying shape whose time 
evolution contains an inherent ambiguity. 

For a simple illustration of this point imagine two coplanar ellipses sharing one focus. 
Suppose they rotate at different rates in their common plane. Let a planet be located at one of 
the intersection points of these ellipses. The values of the elliptic elements needed to describe 
its trajectory would depend upon which ellipse was chosen to parameterise the orbit. Either 
set would be equally legitimate and would faithfully describe the physical trajectory. Thus 
we see that there exists an infinite number of ways of dividing the actual planet's movement 
into motion along its orbit and the simultaneous evolution of the orbit. Although the physical 
trajectory is unique, its description (parametrisation in terms of Kepler's elements) is not. 
A map between two different (though physically-equivalent) sets of orbital elements is a 
symmetry transformation (a gauge transformation, in physicists' jargon). 

Lagrange never dwelled on that point. However, in his treatment (based on direct ap- 
plication of the method of variation of parameters) he passingly introduced a convenient 
mathematical condition which removed the said ambiguity. This condition and possible 
alternatives to it will be the topics of Sections 1 - 3 of this paper. 

1 The real name of the young man invited in 1766 to the Prussian court was Giuseppe Lodovico Lagrangia. 
It was only several years down the road that he became Joseph-Louis Lagrange. 



2 



In 1834 - 1835 Hamilton put forward his theory of canonical transformations. Several 
years later this approach was furthered by Jacobi who brought Hamilton's technique into 
astronomy and, thereby, worked out a new method of deriving the planetary equations 
(Jacobi 1866), a method that was soon accepted as standard. Though the mathematical 
content of the Hamilton- Jacobi theory is impeccably correct, its application to astronomy 
contains a long overlooked aspect that needs attention. This aspect is: where is the Lagrange 
constraint tacitly imposed, and what happens if we impose a different constraint? This issue 
will be addressed in Section 4 of our paper. 



1.2 Osculating Elements vs Orbital Elements 

We start in the spirit of Lagrange. Before addressing the N-particle case, Lagrange referred 
to the reduced 2-body problem, 

•• u r 

? + ^- = > 

(1) 

T — f planet f sun j fJ> G{jTLpi ane i -\- TTl sun ) , 

whose generic solution, a Keplerian ellipse or a hyperbola, can be expressed, in some fixed 
Cartesian frame, as 

r = f{C x ,...,C % ,t) , ¥= i{C u ...,Cs,t) , (2) 

where 

¥) • m 

/ C=const 

Since the problem $IJ) is constituted by three second-order differential equations, its general 
solution naturally contains six adjustable constants Cj. At this point it is irrelevant which 
particular set of the adjustable parameters is employed. (It may be, for example, a set 
of Lagrange or Delaunay orbital elements or, alternatively, a set of initial values of the 
coordinates and velocities.) 

— *■ 

Following Lagrange (1808, 1809, 1810), we employ / as an ansatz for a solution of the 
N-particle problem, the disturbing force acting at a particle being denoted by AF : 2 

F + 4 - = AF . (4) 

Now the "constants" become time dependent: 

r = /(^(t),...,^),*) , (5) 

whence the velocity 

dr _ d£ df_ dQ _ ^ df_ dQ 

dt ~ dt V dQ dt ~ g V da dt ' U 



2 Our treatment covers disturbing forces AF(r, r, t) that are arbitrary vector- valued functions of position, 
velocity, and time. 
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acquires a new input, Y,(®f / dCi)(dCi/ dt) . 

Substitution of (0) and © into the perturbed equation of motion (jlj) leads to three inde- 
pendent scalar second-order differential equations which contain one independent parameter, 
time, and six functions Ci(t) . These are to be found from the said three equations, and 
this makes the essence of the variation-of-parameters (VOP) method. However, the latter 
task cannot be accomplished in a unique way because the number of variables exceeds, by 
three, the number of equations. Thence, though the physical trajectory (comprised by the 
locus of points in the Cartesian frame and by the values of velocity at each of these points) 
is unique, its parametrisation through the orbital elements is ambiguous. This circumstance 
was appreciated by Lagrange who amended the system, by hand, with three independent 
conditions, 

_ df_ dQ 

r dt ' 1 ' 

and thus made it solvable. His choice of constraints was motivated by both physical consider- 
ations and mathematical expediency. Since, physically, the set of functions (Ci(t), ...,Ce(t)) 
can be interpreted as parameters of an instantaneous ellipse, in a bound-orbit case, or an 
instantaneous hyperbola, in a fly-by situation, Lagrange found it natural to make the instan- 
taneous orbital elements Cj osculating. His constraint (J7J) fixes the instantaneous ellipse 
(or hyperbola) defined by (Ci(£), ...,Ce(t)) such that, at each moment of time, it coincides 
with the unperturbed (two-body) orbit that the body would follow if the disturbances were 
to cease instantaneously. This way, Lagrange restricted his use of the orbital elements to 
elements that osculate in the reference frame wherein ansatz (jSJ is employed. Lagrange never 
explored alternative options; he simply imposed (J2J) and used it to derive his famous system 
of equations for orbital elements. 

Such a restriction, though physically motivated, is completely arbitrary from the math- 
ematical point of view. While the imposition of (J7J) considerably simplifies the subsequent 
calculations it in no way influences the shape of the physical trajectory and the rate of 
motion along it. A choice of any other supplementary constraint 

— * 

<fr being an arbitrary function of time and parameters Cj, would lead to the same physical 
orbit and the same velocities. 3 Substitution of the Lagrange constraint (J7J) by its general- 
isation (|SJ) would not influence the motion of the body but would alter its mathematical 
description (i.e., would entail different solutions for the orbital elements). Such invariance 
of a physical picture under a change of parametrisation is called gauge freedom or gauge 
symmetry. It parallels a similar phenomenon well known in electrodynamics and, therefore, 
has similar consequences. On the one hand, a good choice of gauge often simplifies solution 
of the equations of motion. (In Section 3 we provide a specific application to motion in an 
accelerated coordinate system.) On the other hand, a so-called "gauge shift" will occur in 
the course of numerical orbit computation (Murison & Efroimsky 2003). 

Derivation of the conventional Lagrange and Delaunay planetary equations by the VOP 
method incorporates the Lagrange constraint (Brouwer & Clemence 1961). Both systems 

3 In principle, one may endow <& also with dependence upon the parameters' time derivatives of all orders. 
This would yield higher-than-first-order time derivatives of the Cj in subsequent developments requiring 
additional initial conditions, beyond those on f and f, to be specified in order to close the system. We avoid 
this unnecessary complication by restricting €> to be a function of time and the Cj . 
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of equations get altered under a different gauge choice as we now show. The essence of a 
derivation suitable for a general choice of gauge starts with © from which the formula for 
the acceleration follows: 

d 2 r = dj y< dj dCj d$ d 2 f ^ dg_ dQ d4 

dt 2 ~ dt ^ dQ dt dt ~ d 2 t ^ dQ dt dt ' U 

Together with the equation of motion (jlj), it yields: 

W + ~ 2 ~ + E dQ -df + I! = AF ' r=W = \f\ ■ ( 10 ) 

The vector function / was from the beginning introduced as a Keplerian solution to the 
two-body problem; hence it must obey the unperturbed equation ((TJ) . On these grounds the 
above formula simplifies to: 

V & = AF - ^ . (11) 

v da dt dt 1 ' 

This equation describes the perturbed motion in terms of the orbital elements. Together 
with constraint (jHj) it constitutes a well-defined system which may be solved with respect 
to dCi/dt . An easy way to do this is to use the elegant trick suggested by Lagrange: to 

— * 

multiply the equation of motion by df jdC n and to multiply the constraint by — dg/dC n . 
The former operation results in 

(rff dg_dCj\ dl A ^_dfd$ 

8C n \y dCj dt ) 8C n dC n dt ' 1 ' 

while the latter gives 

oc n \ff acj dt ) dc„ y ' 

Having summed these two equalities, we arrive to: 

^ [ ° n ° jJ dt dc n dc n dt oc n ' 1 j 

[C n Cj] standing for the unperturbed (i.e., defined as in the two-body case) Lagrange 
brackets: 

— *■ — p 

df <9g df <9g 

[Cn Cj] = dc~ n dc 3 ~ dc 3 dc~ n • (15) 

— * 

It was agreed above that $ is a function of time and of the "constants" C, , but not 
of their time derivatives. Under this convention, the above equation may be shaped into a 
more convenient form: 
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This is the most general form of the gauge-invariant perturbation equations of celestial 

— * 

mechanics. In the Lagrange gauge, when the $ terms are absent, we can obtain an 
immediate solution for the individual dCi/dt by exploiting the well known expression for 
the Poisson-bracket matrix which is inverse to the Lagrange-bracket one and is offered in 
the literature for the two-body problem. (Be mindful that our brackets (|15|) are defined in 

— * 

the same manner as in the two-body case; they contain only functions / and g that are 
defined in the unperturbed, two-body, setting.) In an arbitrary gauge the presence of the 
term proportional to 8<&/dCj on the left-hand side of complicates the solution for 
dCi/dt, but only to the extent of requiring the resolution of a set of six simultaneous linear 
algebraic equations. 

To draw to a close, we again emphasise that, for fixed interactions and initial conditions, 
all possible (i.e., compatible and sufficient) choices of gauge conditions expressed by the 
vector function <fr lead to a physically equivalent picture. In other words, the real trajectory 
is invariant under reparametrisations permitted by the ambiguity of the choice of gauge. This 
invariance has the following meaning. Suppose the equations of motion for Ci,...,6 , with 
some gauge condition # imposed, render the solution G\ ,...,6(0 • The same equations of 

motion, with another gauge $ enforced, furnish the solution Ci(t) that has a different 
functional form. Despite this difference, both solutions, Cj(t) and Ciif) , when substituted 
back in (jSj), yield the same curve r(t) with the same velocities r(t) . In mathematics this 
situation is called a fibre bundle, and it gives birth to a 1-to-l map of Ci(t) onto Ci(t) , 
which is merely a reparametrisation. In physics this map is called a gauge transformation. 
The entire set of these reparametrisations constitute a group of symmetry known as a gauge 
group. 

— * — * 

Just as in electrodynamics, where the fields E and B stay invariant under gradient trans- 
formations of the 4-potential A 11 , so the invariance of the orbit implements itself through 
the form-invariance of expression (j5J) under the afore mentioned map. The vector r and its 
full time derivative r , play the role of the physical fields E and B , while the Keplerian 
coordinates Ci, C 6 play the role of the 4-potential N 1 . 

A comprehensive discussion of all the above-raised issues can be found in Efroimsky 
(2002, 2003). Interconnection between the internal symmetry and multiple time scales in 
celestial mechanics is touched upon in Newman & Efroimsky (2003). 



2 Planetary equations in an arbitrary gauge 
2.1 Lagrangian and Hamiltonian Perturbations 

We can proceed further by restricting the class of perturbations we consider to those in which 
AF is derivable from a perturbed Lagrangian. This restricted class is still broad enough to 
encompass most applications of celestial-mechanics perturbation theory. 

Let the unperturbed dynamics be described by Lagrangian L (r, r, t) = r /2 — 
U(r , t), canonical momentum p = r and Hamiltonian H (r, p, t) — p 2 /2 + U(r , t) . 
Disturbed motion will be described by the new, perturbed, functions: 

u 2 

L(r, F, t) = L + AL = - U(r, t) + AL(r, ?, t) , (17) 

P = ? + ^ , (18) 
or 
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-> i \ ( d A T, \ ^ 

H(f, p, t) = p? - L = ^- + U + AH , Aif(r, p, i) = - AL - - ( -^-j , (19) 

AH being introduced as a variation of the functional form, i.e., as AH = H(r, p, t) — 
H Q (r, p, t) . The Euler-Lagrange equations written for the perturbed Lagrangian (fT7)) will 
give: 

dU 

f = - + AF , (20) 

where the new term 

d AL d (dAL\ . , 

is the disturbing force. Expression ()21|) reveals that whenever the Lagrangian perturbation 
is velocity independent, it plays the role of the disturbing function. Generally, though, the 
disturbing force is not equal to the gradient of AL but has an extra term generated by the 
velocity dependence. 

Examples in which a velocity dependent AL has been used in a celestial mechanics 
setting include: the treatment of inertial forces in a coordinate system tied to the spin axis 
of a precessing planet (Goldreich 1965); the velocity dependent corrections to Newton's law 
of gravity in the relativistic two-body problem (Brumberg 1992). 



2.2 Gauge- invariant planetary equations 

When the expression (|21j) for the most general force emerging within the Lagrangian for- 
malism is substituted into the gauge- invariant perturbation equation (JHj), it yields: 



^2 \Pn Cj 



df d AL 



d f d ($ + 9Ai ^ 



dt dC n dr dC n dt \~ ' d'v ) dC u 
Since, for a velocity-dependent disturbance, the chain rule 

d AL d AL df d AL d? d AL df dAL d(g + $ 



dC„ 



dr dC„ 



d? dC n 



dr dC„ 



dr 



dC n 



takes place, formula (j22|) may be re-shaped into 



(22) 



(23) 



E [0. o A f 



dAL dAL d$> 



df d 



dt 



dC n 



dr dC n 



dC n dt dC n 



$ + 



d AL \ 
dr J 



(24) 



Next we group terms so that the gauge function $ everywhere appears added to d(AL) / dr 
and bring the term proportional to dCj/dt to the left hand side of the equation: 



E [0. Cj 



df d (dAL 



dC n dCj \ dr 



+ $ 



d 



dC„ 



AL 



1 fdAL\ : 



dr J 



dCj 
~dT 



df d 



dAL d 



dC„ 



dC n dt 



dr dC n 



$ + 



dAV 
dr . 



(25) 



These modifications help us to recognise the special nature of the gauge <fr = — d(AL)/dr 
which will be the subject of discussion in the next subsection. 

Contrast (|25|) with ()16|) : while (|16|) expresses the VOP method in the most general form 
it can have in terms of the disturbing force AF(r, r, t) , equation ()25|) renders the most 
general form in the language of a Lagrangian perturbation AL(r, r, t). 

The applicability of so generalised planetary equations in analytical calculations is com- 
plicated by the nontrivial nature of the left-hand sides of (16) and (25). Nevertheless, the 
structure of these left-hand sides leaves room for analytical simplification in particular situ- 
ations. One such situation is when the gauge is chosen to be 

dAL , 
* = - —^r + t) , (26) 
or 

fj(r, t) being an arbitrary vector function linear in r . (It may be, for example, proportional 
to r or may be equal, say, to a cross product of r by some time-dependent vector.) Under 
these circumstances the left-hand side in (25) reduces to the Lagrange brackets. The situation 
becomes especially simple when dAL/dr happens to be linear in r , in which case we may 
put ff(r, t) = dAL/dr and, thus, employ the trivial Lagrange gauge $ = instead of 
the generalised Lagrange gauge. We shall encounter one such example in section 3.4. 

As already stressed above, the Lagrange brackets are gauge-invariant because functions 
/ and g used in are defined within the unperturbed, two-body, problem (P-EJ) that lacks 
gauge freedom. For this reason, one may exploit, to solve (J2"5)l . the well-known expression 
for the inverse of matrix [Cj Cj] . Its elements are simplest (and are either zero or unity) 
when one chooses as the "constants" the Delaunay set of orbital variables (Plummer 1918): 

d = {L, G, H, M , u, n} 

(27) 



L = y/JTa , G = yfia (1 — e 2 ) , H = y /x a (1 — e 2 ) cos i , 

where /i = G(m sun + m p i anet ) and (e, a, M Q , u, Q , i) are the Keplerian elements: e 
and a are the eccentricity and major semiaxis, M is the mean anomaly at epoch, and the 
Euler angles u> , ft , i are the the argument of pericentre, the longitude of the ascending 
node, and the inclination, respectively. 

The simple forms of the Lagrange and Poisson brackets in Delaunay elements is the proof 
of these elements' canonicity in the unperturbed, two-body, problem: the Delaunay elements 
give birth to three canonical pairs (Qi, Pi) corresponding to a vanishing Hamiltonian: 
(L , — M ) , (G , —uj), (H , — fl) . In a perturbed setting, when only a position-dependent 
disturbing function R(r, t) is "turned on", it can be expressed through the Lagrangian and 
Hamiltonian perturbations in a simple manner, R(r, t) = AL(r, t) = — AH(r, t) , as can 
be seen from the formulae presented in the previous subsection. Under these circumstances, 
the Delaunay elements remain canonical, provided the Lagrange gauge is imposed (Brouwer 
& Clemence 1961). This long known fact can also be derived from our equation (J25)) : if we 
put $ = and assume AL velocity-independent, we arrive to 

™§ = 9 ~M • < 28) 



where 

AL = AL(f(C,t), t) = R(f{C,t), t) = -AH(f{C,t), t) . (29) 
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This, in its turn, results in the well known Lagrange system of planetary equations, provided 
the parameters Cj are chosen as the Kepler elements. In case they are chosen as the 
Delaunay elements, then (|2*5|) leads to the standard Delaunay equations, i.e., to a symplectic 
system wherein the pairs (L , — M Q ) , (G , —uj), (H , — O) again play the role of canonical 
variables, but the Hamiltonian, in distinction to the unperturbed case, no longer vanishes, 
instead being equal to AH = — AL . 

In our more general case, the perturbation depends also upon velocities (and, therefore, 

— * 

AL is no longer equal to — AH ) . Beside this, the gauge $ is set arbitrary. As 
demonstrated in Efroimsky (2002, 2003), under these circumstances the gauge- invariant 
Delaunay-type system is no longer symplectic. However, it turns out that this system regains 
the canonical form in one special gauge, one that coincides with the Lagrange gauge when the 
perturbation bears no velocity dependence. The issue is explained at length in our previous 
papers (Efroimsky & Goldreich 2003a,b). Here we offer a brief synopsis of this study. 

2.3 Generalised Lagrange gauge wherein the Delaunay-type 
system is canonical 

Equation (}2*2*|) was cast in the shape of (J23|) not only to demonstrate the special nature of 
the gauge 

*--^. (3") 

or 

but also to single out the terms in the square brackets on the right-hand side of (|25|): together, 
these terms give exactly the Hamiltonian perturbation. Thus we come to the conclusion that 
in the special gauge (j3T?j) our equation (J23j) simplifies to 

As emphasised in the preceding subsection, the gauge invariance of definition (fT3j) enables us 
to use the standard (Lagrange-gauge) expressions for [C n Cj]^ 1 to get the planetary equations 
from (jSH). 

Comparing (J3TJ) with (|28p. we see that in the general case of an arbitrary AL(r, r, t) one 
arrives from (|3*T|) to the same equations as from ()28|). except that now they contain — AH 
instead of AL. When the orbit is parametrised by the Delaunay variables, those equations 
take the form: 

AT a \ U At T\/T \ a A U 

(32) 
(33) 



dL 


dAH 


d{ - M ) 


dAH 


~dt ' 


d(-M ) ' 


dt 


dL 


dG 


dAH 


d(-u) 


dAH 


~dt 


d(-u) 


dt 


dG 


dH 


dAH 


d{-Q) 


dAH 


~~dt 


d{-Q) ' 


dt 


dH 



(34) 

which is a symplectic system. 4 For this reason we name this special gauge the "generalised 
Lagrange gauge". In any different gauge $ the equations for the Delaunay variables 

In this system H stands not for the Hamiltonian but for one of the Delaunay elements. 
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would contain ^-dependent terms and would not be symplectic. (Those gauge-invariant 
equations, for both Lagrange and Delaunay elements are presented in Efroimsky (2002, 2003) 
and Efroimsky & Goldreich (2003a,b).). This analysis proves the following THEOREM : 
Though the gauge-invariant equations for Delaunay elements are, generally, not 
canonical, they become canonical in the "generalised Lagrange gauge". That this 
Theorem is not merely a mathematical coincidence but has deep reasons beneath it will be 
shown in Section 4 where the subject is approached from the Hamilton- Jacobi viewpoint. 

The above Theorem gives one example of the gauge formalism being of use: an appro- 
priate choice of gauge can considerably simplify the planetary equations (in this particular 
case, it makes them canonical). 

According to ()18j) . the momentum can be written as 

u dAL _ -j dAL 
P = r + = g + $ + , (35) 

or or 

which, in the generalised Lagrange gauge (J3T?|) . simply reduces to 

p = g . (36) 

Vector g was introduced back in ((21-121) to denote the functional dependence of the 
unperturbed velocity upon the time and the parameters Ci . In the unperturbed, two-body, 
setting this velocity is equal to the momentum canonically conjugate to the position r (this 
is obvious from (|T%|) . for zero AL ). This way, in the unperturbed case equality (plo]) is 
fulfilled trivially. The fact that it remains valid also under perturbation means that, in the 
said gauge, the canonical momentum in the disturbed setting is the same function of time 
and "constants" as in the unperturbed, two-body, case. Thus we have established that the 
instantaneous Keplerian ellipses (hyperbolae) defined in gauge (|HU|) osculate the trajectory 
in phase space. 

Not surprisingly, the generalised Lagrange gauge (|3Ti|) reduces to $ = in the simple 
case of velocity-independent disturbances. 

3 Gauge Freedom and Freedom of Frame Choice 

3.1 Osculating ellipses described in different frames of reference. 

The essence of the VOP method in celestial mechanics is the following. A generic two-body- 
problem solution expressed by 

r = f(C,t) , (37) 
= i(C,t) , (38) 





1L J- 

f 2 f 



(39) 



is employed as an ansatz to solve the disturbed problem: 

r = f(C(t), t) , (40) 
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df df dd 



- dg dg dd d$> 

r ~ ~dt + dCl ~dt + ~dt 

(42) 

P f dd dt dt ' 

As evident from (JHJ), our choice of a particular gauge is equivalent to decomposition of 
the physical motion into a movement with velocity g along the instantaneous ellipse (or 
hyperbola, in the fly- by case), and a movement associated with the ellipse's (or hyperbola's) 
deformation that goes at the rate <fr . It is then tempting to state that a choice of gauge is 
equivalent to a choice of an instantaneous comoving reference frame wherein to describe the 
motion. Such an interpretation is, however, incomplete. Beside the fact that we decouple 

— *■ 

the physical velocity in a certain proportion between g and <fr , it also matters which 
physical velocity (i.e., velocity relative to what frame) is decoupled in this proportion. In 
other words, our choice of the gauge does not yet exhaust all freedom: we can still choose in 
what frame to write ansatz fl40(l . We may write it in inertial axes or in some accelerated 
system. For example, in the case of a satellite orbiting an accelerated and precessing planet 
it is convenient to write the ansatz for the planet-related position vector. 

The above kinematic formulae (}4T)j) - (pf2"j) do not yet contain information about our choice 
of the reference system in which we implement the VOP method. This information shows 
up at the next stage, when expression (|42j) is inserted into the dynamical equation of motion 
r = — (fxr/r 3 ) + AF to yield 

dg^dQ dS f dAL _ d_ fdAL\ 
DCi dt dt dr dt { d? ) ' 

Complete information about the reference system in which we put the VOP method to work 
(and, therefore, in which we define the orbital elements Cj ) is contained in the expression 
for the perturbation force AF . For example, if the operation is carried out in an inertial 
coordinate system, AF contains physical forces solely. However, if we wish to implement 
the VOP approach in a frame moving with a linear acceleration a , then AF also contains 
the inertial force — a . In case this coordinate system rotates relative to inertial ones at a 
rate jl , then AF also includes the inertial terms — 2/Xx r — p, x r — /xx(/7xr). In 
considering the motion of a satellite orbiting an oblate precessing planet it is most reasonable, 
though not obligatory, to apply the method (i.e., to define the time derivative) in axes that 
precess with the planet. However, this reasonable choice of coordinate system still leaves us 
with the freedom of gauge nomination. 

3.2 Relevant Example 

Gauge freedom of the perturbation equations of celestial mechanics finds an immediate 
practical implementation in the description of test particle motion around an precessing 
oblate planet (Goldreich, 1965). It is trivial to extend this to account for acceleration of the 
planet's centre of mass. 
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Our starting point is the equation of motion in the inertial frame 



dU 
dr 



(44) 



where U is the total gravitational potential and time derivatives in the inertial axes are 
denoted by primes. Suppose that the planet's spin axis precesses at angular rate fi{t) and 
that the acceleration of its centre of mass is given by a(t). In a coordinate system attached 
to the planet's centre of mass and precessing with it, inertial forces modify the equation of 
motion so that it assumes the form: 



dU 
dr 



2/t x r — fix r — fi x (fi x r) 



(45) 



time derivatives in the accelerated frame being denoted by dots. 

To implement the VOP approach in terms of the orbital elements defined in the accel- 
erated frame, we note that the disturbing force on the right-hand side of (|45jl is generated 
according to (J2~T|) by 



AL r, r, t 



R + r(fi x r) 



x rj-(/Lt x rj 



a • r 



(46) 



where we denote by R(r , t) the gravitational-potential perturbation (which the perturbation 
of the overall gravitational potential U). Since 



d AL 
di 



/x x r 



(47) 



the corresponding Hamiltonian perturbation reads: 



AH 



1 dAL 
AL + - — — 

2 I dr 



= — [R + p • (/X x r) — a-r 



[R + (r x p) • p, 



(48) 



— a-r 



r x p being the satellite's orbital angular momentum in the inertial 



with vector J 
frame. 

According to (|35|) and (J47|) . the momentum can be written as 

p = g + # + /xx/ , 
whence the Hamiltonian perturbation becomes 



AH 



R 



f x g /2 



fix f) ■ (fix f) - a • / 



(49) 



(50) 



3.3 Elements defined in an accelerated, rotating frame 
that osculate in the comoving inertial frame 

In this subsection we recall a calculation carried out by Goldreich (1965) and Brumberg, 
Evdokimova & Kochina (1971) and demonstrate that it may be interpreted as an example 
of nontrivial gauge fixing. 
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Let us implement the VOP method in a frame that is accelerating at rate a and rotating 
at angular rate /X relative to some inertial system S. This means that, in the VOP equation 
flMl), AL is given by formula flH| and AH by 1(5(1)1 . 

We now choose to describe the motion in the generalised Lagrange gauge (|3T)j). so the 
expression + ft x r ) on the right-hand side of (|5()jl vanishes (as follows from (|47)l). 

— *■ 

and the expression for AH in terms of / and g has the form: 



AH 



R{f,t)+ix-{fxi) - a / . (51) 



At the same time, the generic expression for the VOP given by (J25|) simplifies to (|31|) . 
Insertion of (|5Tj) therein leads us to 

dC f) 

[CrCi] ~df = da. + /x -^ xg ) - s-/] • (52) 

Interestingly, this equation does not contain p, even though it is valid for non-uniform 
precession. 

As explained in subsection 2.3, in the generalised Lagrange gauge the vector g is equal 
to the canonical momentum p = r + dAL/dr. In the case when the velocity dependence 
of AL is called into being by inertial forces, the momentum is, according to (|47j) . 

u d AL u _ _ 

p = r + — = r + fi x r , (53) 

dr 

which is the particle's velocity relative to the inertial frame comoving with the accelerated, 
rotating frame. In this sense we may say that our elements are defined in the accelerated, 
rotating frame but osculate in the comoving inertial one. 

In the appendix we provide explicit expression for each of the partial derivatives of p, ■ J 
that appears in the planetary equations (|52|). 



3.4 Elements defined in the accelerated, rotating frame, 
that osculate in this frame 

Here we not only define the elements in the accelerated, rotating frame, but we also make 
them osculate in this system, i.e., we make them satisfy $ = 0. In this gauge, expression 
(|50j) takes the following form: 

AH = - [R(f, t) + f2-(f xg) + (jtxf)-(jtxf) - a-/] . (54) 

while equation (|25[). after some algebra, 5 looks like this: 

da dAH 



dt ' oa 



(55) 



5 Due to (|47l) , the second term on the left-hand side in (|25l) is proportional to [d(fl x /)/ dCj] Cj = /7x $ 
and, therefore, vanishes. The second term on the right-hand side simplifies in accordance with the simple 
rule A ■ (B x C) = (C x A) • B . (See Efroimsky & (Goldreich 2003b).) 
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When substituting (J54j) into (J55|) . it is convenient to rent the expression for AH apart and 
to group the term (ft x f) ■ (ft x f) with the last term on the right-hand side of (jo3j) : 



(56) 



In so writing Ij56j) we have deliberately cast it into a form that eases comparison with fl52(l . 

In the appendix we set up an apparatus from which the partial derivatives of the inertial 
terms with respect to the orbital elements may be obtained. We also show that some of these 
derivatives vanish. However, a complete evaluation of the inertial input to the planetary 
equations in the ordinary Lagrange gauge involves a long and tedious calculation which we 
do not carry out. 



3.5 Comparison of the two gauges 

One of the powers of gauge freedom lies in the availability of gauge choices that simplify the 
planetary equations, as we can see from contrasting (}52|) with (}56|) . While the latter equation 
is written under the customary Lagrange constraint (i.e., for elements osculating in the frame 
where they are defined), the former equation is written under a nontrivial constraint called 
the "generalised Lagrange gauge." The simplicity of ()52|) speaks for itself. 

By identifying the parameters Ci with the Delaunay variables, one arrives from (p)2"j) and 
flST)j) to the appropriate Delaunay-type equations (see Appendix I to Efroimsky & Goldreich 
2003a). The Lagrange equations corresponding to (J52|) and to (J56|) may be derived from 
each of these two equations by choosing Cj as the Kepler elements and using the appro- 
priate Lagrange brackets. These Lagrange equations are written down in the Appendix to 
Efroimsky & Goldreich (2003b) to which we refer the interested reader. 

Although the planetary equations are much simpler in the generalised Lagrange gauge 
than in the ordinary Lagrange gauge, some of these differences are less important than 
others. In many physical situations, though not always, the ft 2 and ft terms in (jo^|) are of 
a higher order of smallness compared to those linear in ft , and therefore may be neglected, 
at least for sufficiently short times. 6 



4 Planetary Equations and Gauges in the Hamilton- 
Jacobi Approach 

In this section we demonstrate that the derivation of planetary equations in the N-particle 
(N > 3) case, performed through the medium of Hamilton- Jacobi method, implicitly con- 
tains a gauge-fixing condition not visible to the naked eye. We present a squeezed account of 
our study; a comprehensive description containing technical details may be found in Efroim- 
sky & Goldreich (2003a). 

The Hamilton- Jacobi analysis rests on the availability of different canonical descrip- 
tions of the same physical process. Any two such descriptions, (q, p, H(q,p)) and 

6 As an example of an exception to this rule, we mention Venus whose wobble is considerable. This means 
that, for example, the /2 term cannot be neglected in computations of circumvenusian orbits. 
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(Q, P, H*(Q, P)) , correspond to different parametrisations of the same phase flow, and 
both obey Hamilton's equations. Due to the latter circumstance the infinitesimally small 
variations 

dO = pdq - Hdt (57) 

and 

d9 = PdQ - H*dt (58) 
are perfect differentials, and so is their difference 

- dW = d6 - d6 = PdQ - pdq - (H* - H) dt . (59) 

Here, vectors q, p, Q, and P each contain N components. Given a phase flow parametrised 
by a set (g, p, H(q,p,t)) , it is always useful to simplify the description by a canonical 
transformation to a new set (Q, P, H*(Q, P,t) ) , with the new Hamiltonian H* being 
constant in time. Most advantageous are transformations that nullify the new Hamiltonian 
H* , because then the new canonical equations render the variables (Q, P) constant. 
A powerful method of generating such transformations stems from (J59|) being a perfect 
differential. It is sufficient to consider W to be a function of the time and only two other 
canonical variables, for example q and Q . Then (|59|) may be written as 

dW dW dW 

— —— dt — — — dQ — dq = PdQ - pdq + (H - H*) dt (60) 

at oQ oq 

from which it follows that 

dW dW dW 

P = -9Q • "=17 ' H {q ,„, t) + ljr ^H- ( Q,P, t) . (61) 

Inserting the second equation into the third and assuming that H*(Q, P, t) is simply a 
constant, we get the famous Jacobi equation 

/ dW \ dW 

whose solution furnishes the transformation-generating function W. The elegant power of the 
method becomes most visible if the constant H* is set to zero. Under this assumption the re- 
duced two-body problem is easily resolved. Starting with the three spherical coordinates and 
their canonical momenta as (q, p) , one arrives to canonically conjugate constants (Q, P) 
that coincide with the Delaunay elements (|27jl: (Qi, Pi) = (L , — M a ) ; (Q2, P2) = 
(G, -u) ; (Q 3 , P 3 ) = (H, -fi) . 

Extension of this approach to the N-particle problem begins with consideration of a dis- 
turbed 2-body setting. The number of degrees of freedom is still the same (three coordinates 
q and three conjugate momenta p), but the initial Hamiltonian is perturbed: 

. = d(H + AH) = _ d(H + AH) 

dp ' dq 

While in f!60|) - one begins with the initial Hamiltonian H and ends up with H* = , 
the method may be extended to the perturbed setting by accepting that now we start 
with a disturbed initial Hamiltonian H + AH and arrive, through the same canonical 
transformation, to an equally disturbed eventual Hamiltonian H* + AH = AH . Plugging 
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these new Hamiltonians into (|6(Jj) leads to cancellation of the disturbance AH on the 
right-hand side, whereafter one arrives to the same equation for W(q, Q, t) as in the 
unperturbed case. Now, however, the new canonical variables are no longer conserved but 
obey the dynamical equations 

9A£ _ dAH 

V dP • dQ ' [ ' 

Because the same generating function is used in the perturbed and unperturbed cases, the 
new, perturbed, solution (q, p) is expressed through the perturbed "constants" Q(t) and 
P(t) in the same manner as the old, undisturbed, q and p were expressed through the old 
constants Q and P . This form-invariance provides the key to the N-particle problem: one 
should choose the transformation-generating function W to be additive over the particles 
and repeat this procedure for each of the bodies, separately. 

Armed with this preparation, we can proceed to uncover the implicit gauge choice made 
in using the Hamilton- Jacobi method to derive evolution equations for the orbital elements. 
To do this we substitute the equalities 

• dAH dAH dq_ dAH dp_ 

Q ~ ~dP~ ~ ~dq~ dP + ~dp~ dP ^ ' 



and 

p dAH dAH dq dAH dp 

dQ dq dQ dp dP 

into the expression for the velocity 



(66) 



dq dq ■ dq ■ 
q= d-t + dQ Q+ dP P - (67) 



This leads to 



dq ( dq dq dq dq\ d AH ( dq dp dq dp\ d AH 
q = di + [dQ dP ~ dP dQ J dq + [dQ dP ~ dP dQ ) dp 

(68) 

fdAH\ dq 

= o + W) q:t • o^m- 

where we have taken into account that the Jacobian of the canonical transformation is unity: 

dq_ dp_ dq_ dp_ 

dQ dP dP dQ 1 ' 

To establish the link between the regular VOP method and the canonical treatment, compare 
(JOB} with (jHJ). We see that the symplectic description necessarily imposes a particular gauge 
$ = dAH I dp. 

It can be easily demonstrated that this special gauge coincides with the generalised 
Lagrange gauge (J3Uj) discussed in subsection 2.2. To that end one has to compare the 
Hamilton equation for the perturbed Hamiltonian ((EH) , 

d (H + AH) dAH , . 

dp dp 
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with the definition of momentum from the Lagrangian (|17|) . 



d(L(q, q, t) + AL(q, q, t)) . d AL 

p = M = q + ~w ■ (71) 

Equating the above two expressions immediately yields 

fdAH\ (dAL\ , . 

which coincides with (J3U|) . Thus the transformation generated by W(q, Q, t) is canonical 
only if the physical velocity q is split in a fashion prescribed by (|5Sjl. i.e., if (J72*Jl is fulfilled. 
This is exactly what our Theorem from subsection 2.2 says. 

To summarise, the generalised Lagrange constraint, $ = — d AL/dq , is tacitly in- 
stilled into the Hamilton- Jacobi method. Simply by employing this method (at least, in 
its straightforward form), we automatically fix the gauge. 7 By sticking to the Hamiltonian 
description we sacrifice gauge freedom. 

Above, in subsection 2.3, we established that in the generalised Lagrange gauge the 
momentum coincides with g . We now can get to the same conclusion from ([68)1. (fTTj) and 
(E2): 

s d(L(q, q,t) + AL(q,q,t)) =q _^ = g . (?3) 
oq 

Thus, implementation of the Hamilton- Jacobi theory in celestial mechanics demands the 
orbital elements to osculate in phase space. Naturally, this demand reduces to that of 
regular osculation in the simple case of velocity-independent AL . 



5 Conclusions 

In the article thus far we have studied the topic recently raised in the literature: the planetary 
equations' internal symmetry that stems from the freedom of supplementary condition's 
choice. The necessity of making such a choice constrains the trajectory to a 9-dimensional 
submanifold of the 12-dimensional space spanned by the orbital elements and their time 
derivatives. Similarly to the field theory, the choice of the constraint (= the choice of gauge) 
is vastly ambiguous and reveals a hidden symmetry instilled in the description of the N-body 
problem in the language of orbital elements. 

We addressed the issue of internal freedom in a sufficiently general setting where a per- 
turbation to the two-body problem depends not only upon positions but also upon velocities. 
Such situations emerge when relativistic corrections to Newton's law are taken into account 
or when the VOP method is employed in rotating systems of reference. 

Just as a choice of an appropriate gauge simplifies solution of the equations of motion 
in electrodynamics, an alternative (to that of Lagrange) choice of gauge in the celestial me- 
chanics can simplify orbit calculations. We provided one such example, a satellite orbiting a 
precessing planet. In this example, the choice of the generalised Lagrange gauge considerably 
simplifies matters. 

7 An explanation of this phenomenon from a different viewpoint is offered in Section 6 of Efroimsky (2003) 
where the Delaunay equations are derived also through a direct change of variables. It turns out that the 
outcome retains the symplectic form only if an extra constraint is imposed by hand. 
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We have explained where the Lagrange constraint tacitly enters the Hamilton- Jacobi 
derivation of the Delaunay equations. This constraint turns out to be an inseparable (though 
not easily visible) part of the method: in the case of momentum-independent disturbances, 
the N-body generalisation of the 2-body Hamilton- Jacobi technique is legitimate only if we 
use orbital elements that are osculating. In the situation where the disturbance depends 
not only upon positions but also upon velocities, another constraint (which we call the 
"generalised Lagrange constraint") turns out to be stiffly embedded into the Hamilton- Jacobi 
development of the problem. 

Unless a specific constraint (gauge) is imposed by hand, the planetary equations assume 
their general, gauge- invariant, form. In the case of a velocity-independent disturbance, any 
gauge different from that of Lagrange drives the Delaunay system away from its symplectic 
form. If we permit the disturbing force to depend also upon velocities, the Delaunay equa- 
tions retain their canonicity only in the generalised Lagrange gauge. Interestingly, in this 
special gauge the instantaneous ellipses (hyperbolae) osculate in phase space. 

Briefly speaking, N-body celestial mechanics, expressed in terms of orbital elements, is 
a gauge theory but it is not strictly canonical. It becomes canonical in the generalised La- 
grange gauge. 
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Appendix 

In this appendix we set up an apparatus from which one may evaluate the partial deriva- 
tives with respect to the orbital elements of inertial terms that appear in the planetary 
equations derived in sections 3.3 and 3.4. We then show that some of these derivatives van- 
ish. Following that, we derive explicit expressions for each derivative of p, • x gj , which 
provides a complete analytic evaluation of the rotational input in the generalised Lagrange 
gauge. The topic is further developed (and the appropriate generalised Lagrange system of 
equations is presented) in (Efroimsky & Goldreich 2003b). 

To find the explicit form of the dependence / = /(Cj, t) , it is conventional to introduce 
an auxiliary set of Cartesian coordinates q , with an origin at the gravitating centre, and 
with the first two axes located in the plane of orbit. The q coordinates are easy to express 
through the major semiaxis a , the eccentricity e and the eccentric anomaly E : 

q ± = a(cosE — e) , q 2 = aVl — e 2 sinE , g 3 = , (74) 

where E itself is a function of the major semiaxis a , the eccentricity e , the mean anomaly 
at epoch, Mq , and the time, t . The time dependence is realised through the Kepler equation 

E - e sin£ = M , (75) 

where 

M = M a + /i 1/2 / a- 3/2 dt . (76) 



18 



The inertial-frame-related position of the body reads: 

r — f (ft, i, lj, a, e, M Q ; t) = R(fl, i, uo) q(a, e, _E(a, e, M D , t) ) , (77) 

R(fi, i, cj) being the matrix of rotation from the orbital-plane-related coordinate system 
q to the fiducial frame (x, y, z) in which the vector r is defined. This rotation is 
parametrised by the three Euler angles: inclination, i ; the longitude of the node, Q ; and 
the argument of the pericentre, u . 

In the unperturbed two-body setting the velocity is expressed by 

g = ^ f(Q, i, u, a, e, M ; t) = ( ^ J R(fi, i, u) ( || J . (78) 

V / a, e, M V / a, e 

— * 

One can similarly calculate partial derivatives of / with respect to M a : 

° f(Q, i, u } a, e, M o] t)= ) R(Q, i, u) ( ^) , (79) 



0M o " v ' ' ' ' ' u ' ' \dMj . _ v ' ' ' 



a, e, £ 



— » 

whence it becomes evident that df/dM Q is parallel to g and, hence, 



9M, 







fi, i, oj, a, e, t 



By a similar trick it is possible to demonstrate that d(f x g)/dM Q is proportional to 

— *■ — * 

d(f x gj/dE and, therefore, to <9(/ x g)/dt . Hence, this derivative vanishes (because in 

— * 

the two-particle case the cross product / xg is an integral of motion). This vanishing of 
9(f x g)/dM a , along with (jS0|). implies: 

r x ( H) = ° ■ (8i) 

In the situation when the parameters Cj are implemented by the Delaunay elements, a 
similar sequence of calculations leads to 



df_ 

dM Q 



= (82) 

Q, u>, L, G, H, t 



and, appropriately, to 



/' x ( M] = • (83) 

n,cv,L,G,H,t 



0M o 



We can proceed much farther in the generalised Lagrange gauge, at least in so far as 
derivatives of the rotational input J/fi = fxg are concerned. (We remind that here and 
everywhere \i stands for the reduced mass, while p, denotes the precession rate.) 

As we proved above, this cross product is independent of M a and, hence, 
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Since J is orthogonal to the orbit plane, it is invariant under rotations of the orbit 
within its plane, whence 

-0. ^ (85, 

To continue, we note that in the two-body setting the ratio 3/fi , is known to be equal to 



a(l — e 2 ) w where w is a unit vector perpendicular to the unperturbed orbit's plane. 
Moreover, in planet-associated noninertial axes (x, y, z) with corresponding unit vectors 
(x, y, z) , the normal to the orbit is expressed by 



w = x sin % sin Q — y sin i cos Q + z cos % . (86) 



Hence, 



and 



. d (fx g) _ , d (ygj - e 2 ) ) i rr^ 



_ d (fx g) _ ^ a (ygj - e 2 ) ) 

" ' ar^ = " ' w — - = " V ^ ' (88) 



where fx± = /j, x sin i sin Q — fi y sin z cos f2 + fi z cos i is the orthogonal-to-orbit 
component of the precession rate. The remaining two derivatives look: 



d (f x g) i dw 



(89) 



and 



\J a (1 — e 2 ) {yU x sin z cos f2 + /x y sin % sin £7} 



^ g (fx g) r— _ dw 

M • -^-^ = y/a (1 - e 2 ) M • — = 

(1 — e 2 ) {/^ cos? sinf2 — /i^ cos? cosf2 — jjL z sin «'} 



(90) 



As for the derivatives of a • / , they may be calculated directly from the expression for 
f(Q, cu, i, a, e, M a ; t) presented above. However, the resulting expressions are cumbersome 
so we do not present them here. 
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